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ABSTRACT 

We calculate the location and stability of the and L5 Lagrange equilibrium points in the circular 
restricted three-body problem as the binary system evolves via gravitational radiation losses. Relative 
to the purely Newtonian case, we find that the L4 equilibrium point moves towards the secondary 
mass and becomes slightly less stable, while the L5 point moves away from the secondary and gains 
in stability. We discuss a number of astrophysical applications of these results, in particular as a 
mechanism for producing electromagnetic counterparts to gravitational-wave signals. 
Subject headings: black hole physics - relativity - gravitational waves - galaxies: nuclei 
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1. INTRODUCTION 

Recent advances in numerical relativity (NR) have led 
to an increasing interest in the astrophysical implica- 
tions of black hole (BH) mergers. Of particular interest 
is the possibility of a distinct, luminous electromagnetic 
(EM) counterpart to a gravitational wave (GW) signal 
(IBloom et al.l 120091 ). If such an EM counterpart could 



be identified with a LISAQ detection of a supermassive 
BH binary in the merging process, then the host galaxy 
could likely be determined. A large variety of potential 
EM signatures have recently been proposed, almost all 
of which require some significant amount of gas in the 
near vicinity of the merging BHs. 

Gas in the form of accretion disks around single mas- 
sive BHs is known to produce some of the most luminous 
objects in the universe. However, very little is known 
about the behavior of accretion disks around two BHs, 
particularly at late times in their inspiral evolution. In 
Newtonian disks, it is believed that a circumbinary accre- 
tion disk will have a central gap of much lower density, ei- 
ther preventing acc retion altogether, or at least decreas- 
ing it signif icantly (jPringld [l991l : lArtvmowicz fc Lubowl 
11994 119961 ). When including the evolution of the bi- 
nary due to GW losses, the BHs may also decouple 
from the disk at the point when the GW inspiral time 
becomes shorter than the gaseous inflow time at th e 
inner edge of the disk (jMilosavlievic fc Phinnevll2005l) . 
This decoupling should effectively stop accretion onto 
the central object until the gap can be filled on an inflow 
timescale. However, other semi-analytic calculations pre- 
dict an enhancement of accretion power as the evolving 
binary squeezes the gas around the primary BH, leading 
to a rapid increase in luminosity shortly before merger 
(jArmitage fc Nataraianll2002t iChang et al.ll2009f ). 

Setting aside for the moment the question of how the 
gas can or cannot reach the central BH region, a num- 
ber of recent papers have shown that if there is suf- 
ficient gas present, then an observable EM signal is 
likely. iKrolikl (|2010[ ) used analytic arguments to esti- 
mate a peak luminosity comparable to that of the Ed- 
dington limit, independent of the detailed mechanisms 
for shocking and heating the gas. Using relativistic 

^ |http://lisa. nasa.gov I 



magneto-hydrodynamic simulations in 2D. lO'Neill et al.) 
(200!|) showed that the prompt mass loss due to GWs 
actually leads to a sudden decrease in luminosity fol- 
lowing the merger, as the gas in the inner disk tem- 
porarily has too much energy and angular momen- 
tum to accrete efficiently. Full NR simulations of the 
final few orbits of a merging BH binary have now 
been carried out including th e presenc e of EM fields 
m a vacuum (iPalenzuela et al.l [2009: Mo sta et al.l 120101 : 
iPalenzuela et al."2010) an d also gas, treated as test par- 
ticles in .van Meter et a l ] (I2010D and as an ideal fluid in 
iBode et al.l ()2010D and iFarris et al.l (|2010). The simula- 
tions including matter all suggest that the gas can get 
shocked and heated to high temperatures, thus leading 
to bright counterparts in the event that sufficient gas is 
in fact present. One estimate of "sufficient" is the con- 
dition that the gas be optically thick to electron scat- 
tering, which will lead to Eddington-rate luminosities 
(IKrolikll2010l) . Thus we require npR > 1, where k — QA 
cm^g~^ is the opacity and p is the characteristic density 
over a length scale R ^ IOGA/bh/c^- Then the total 
mass is a modest Mgas ~ pR^ w 5 x 10^^ g, where 
M^ = Mbh/(1O'^M0). More gas will generally lead to 
longer, but not necessarily more luminous signals. 

The motivation for this work is to understand poten- 
tial dynamical processes that could lead to a significant 
amount of gas around the BHs at the time of merger. We 
focus on stable orbits in the circular, restricted three- 
body problem, i.e., a system with two massive objects 
orbiting each other on circular, planar Keplerian orbits, 
with the third body being a test particle free to move 
out of the plane. Gas or stars that get captured into 
these stable regions at large binary separation may re- 
main trapped as the binary orbit shrinks due to radia- 
tion reaction. In this paper, we integrate the restricted 
three-body equations of motion in the corotating frame 
for a range of BH mass ratios, including the evolution 
of the binary with leading-order post-Newtonian GW 
losses. With the exception of these non-conservative cor- 
rections, all dynamics are considered at the purely New- 
tonian level. 

We find that particles initially close to the £4 and L5 
Lagrange points remain on stable orbits throughout the 
adiabatic evolution of the BH binary. This result imme- 
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diately leads to a number of astrophysical predictions. 
If diffuse gas is captured from the surrounding medium, 
and irradiated by an accretion disk around either BH, 
it may produce strong, relatively narrow emission lines 
with large velocity shifts, varying on an orbital time 
scale. If individual stars are captured, they may even- 
tually get tidally disrupted by one of the BHs shortly 
before merger, providing an enormous source of material 
for accretion. Stars may also get ejected from the system 
at late times, producing ultra-hyper-velocity stars mov- 
ing away from the galactic center at a significant fraction 
of the speed of light. Clusters of stars could also get cap- 
tured or form in situ at large binary separation, then get 
compressed adiabatically as the BH system shrinks, ulti- 
mately collapsing from gravitational instabilities. Com- 
pact objects such as neutron stars or stellar-mass BHs 
would not get tidally disrupted in most cases, but could 
still possibly be detected with LISA as perturbations to 
the GW signal. 

While this manuscript w as in preparation, a similar 
work bv lSeto fc Miltol (I2OIOI) was published, with some of 
the same conclusions. One major distinction is that they 
included first-order post-Newtonian terms in their equa- 
tions of motion, while we do not. Also, they "have not 
found notable qualitative differences" between L4 and 
L5, while we do. Many of our astronomical predictions 
are qualitatively similar, but have been derived indepen- 
dently. 

The outline for this paper is as follows: in Section [2] 
we present the three-body equations of motion, including 
GW losses, and analyze the stability and location of the 
equilibrium Lagrange points. These analytic predictions 
are tested in Section [3] with numerical simulations of a 
large number of test particles. In Section|4]we explore in 
more detail some potential astronomical applications of 
these results, and in Section [5] we conclude. 

2. LOCATION AND EVOLUTION OF STABLE LAGRANGE 

POINTS 

We consider a three-body system wherein the larger 
bodies have masses Mi and M2 (Mi > M2; M = 
Ml + M2) and the third object is a test particle of 
negligible mass niQ. We define dimensionless masses 
Hi = Mi/M and ^2 = M2/M. Mi and M2 move on 
circular, Keplerian orbits around the center of mass at 
the origin with semi- major axis a{t). Unless otherwise 
stated, we use geometrized units with G — M = c = 1. 
In the corotating frame we define coordinates {x,y,z), 
with the primary BH located at (— /i2 a, 0, 0) and the sec- 
ondary at (/ii a, 0, 0). The angular velocity of the system 
in the lab frame is n = a~^^'^. The binary is orbiting 
in the direction for the standard definition of spher- 
ical coordinates. The equations of motion for the test 
particle may be written as 



X = 2ny + n x + hy 
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Figure 1. Contours of the Jacobi constant Cjo for test particles 
at rest in the rotating frame, for secondary mass = 0.02. The 
stable Lagrange points L4 and L5 are located at local minima of 
Cjo, and form equilateral triangles with the two massive bodies. 
The co-linear Lagrange points Li, L2, and L3 are located at saddle- 
points of Cjo and are unstable equilibrium points. The x— and 
y— axes have been normalized to the binary separation a. 



where rj and r2 are the distances to the primary and sec- 
ondary bodies, respectively. These reproduce the classi- 
cal Newtonian system in the limit of ri — 0, i.e. when 
gravitational radiation reaction is turned off. 

It is well-known from classical mechanics that there 
exist five equilibrium points where a test particle can re- 
main at r est in the corota ting frame (x = j/ = z = 0) 
([Murray fc DermottI 119990 . These points are known as 
the Lagrange equilibrium points {L\ . . . L5), and are plot- 
ted in Figure [1] for a binary with mass ratio fi2 = 0.02. 
The CO- linear points Lj, L2, and L3 are always unstable, 
while the triangular p oints L4 and Lr are hnear ly stable 
only for fi2 < 0.0385 (|Murrav fc DermottI [T999h . In the 
Newtonian problem, there is a single integral of motion, 
the Jacobi constant: 



^ _ 2/ 2 I 2\ , o ( ^1 I 

C] =n (x +y ) + 2[ \ 



"2 »-.2 ~2 (2) 



X — y — z . 



Setting the particle's velocity to zero in the corotating 
frame, the contours of constant Cjo = Cj{v = 0) in the 
X — y plane define "zero velocity curves," as plotted in 
Figure [TJ L4 and L5 are located at local minima of Cjo, 
while Li, L2, and L3 are saddle-points (note that Cj has 
the form of a negative energy, so that local minima are 
actually points of maximum potential energy, yet are still 
stable to small perturbations because of the restoring 
Coriolis forces in the rotating frame). One well-known 
implication of this feature is that gas bound to cither 
All or M2 can pass through Li and subsequently accrete 
onto the other object, a process known as mass trans- 
fer in standard stellar binary evolution. Additionally, 
any particle with Cj > Cjq^L^) would be excluded from 
crossing into the tadpole-shaped regions around L4 and 



Murray! (|1994f) explore the effects of generalized non- 
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gravitational forces acting on test particles in the re- 
stricted three-body problem, including a detailed study 
of how the position and stability of the Lagrange equi- 
librium points are affected. At first glance, our problem 
woul d se em quite similar, with the n terms in equations 
(fTa |) . (frbl) acting as psuedo-drag forces: as the binary or- 
bit shrinks and n increases, test particles would feel like 
they are accelerating in the —(j) direction in the coro- 
tating frame , thus slowing down in the inertial frame. 
In this case, iMurrayl (|1994l ) show that L4 and L5 both 
move in the -|-(/) direction, thus compensating for the in- 
ertial drag and still maintaining x — y — 0. However, 
when consideri ng t he time-varying nature of n and a in 
equations p^ - (jldp due to GW losses, we find that this 
drag force analog y br eaks down. Numerical integration 
of equations (fTa|) - (jldp show that in fact, the opposite oc- 
curs, with L4 and L5 moving in the —(j) direction in the 
rotating frame. 

Remarkably, this result can still be understood via an 
inertial drag term, if one applies it not to the test par- 
ticle, but to the secondary mass. In order to study the 
evolution of the Trojan asteroids under the influence of a 
radial migration by Jupiter, Fleming & Hamilton (2000) 
apply an artificial drag force that acts only on Jupiter 
in the form F = — fcv in the inertial frame. Of course, 
this causes the secondary to lose energy and angular mo- 
mentum and thus orbit faster, increasing n as in the BH 
binary case. Since there is no requirement on the rela- 
tive masses of the secondary M2 and the test particle mg 
to ensure stability of L4 and L5 (iSalo fc Yoder 1988i). we 
can in fact adopt the drag force analysis of lMurravr([l994f ) 
to include GW losses by reversing the roles of the sec- 
ondary BH and the test particle. In the rotating frame, 
the inertial drag force is given hy = —k(x — ny) and 
Fy = —k{y + nx). For station ary particles at the equi- 
librium points, equations ()lap . (|lb|) become 



kny = —r?x 
knx = n^y — 



M2- 



a/ii 



^3 



Ai2 



y, 



(3a) 
(3b) 



where k = k/iriQ, and we have replaced the GW acceler- 
ation terms fiy and fix with those appropriate for drag 
forces Fx /mo and Fy/mQ. 

In order to reproduce the GW evolution of equation 
(jldp. we must choose the magnitude of the drag force 
appropriately. This is done by equating the energy loss 
due to the shrinking orbit with that of the inertial drag 
force. For circular Keplerian orbits, we have 



F V 



and 



-1/3^ 



(4) 



(5) 



Equating and ([5]), along with pd)) allows us to solve 
for k: 

- In 32 o/q 
3 n 5 

Following lMurravl (|1994[ ) , we multiply equation (|3a|) by 



y and add to (I3bp times x, giving 

knr'^ = Hiii2ay \ \ 



(7) 



In the limit of small fi2, o- • 
points and we can write 

K ~ — SlIKpL 
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ri at the stable Lagrange 



(8) 



(2-2cOS0l)3/2 



where ipL is the location of the Lagrange point in the 
presence of an inertial drag force. When considering the 
orbital evolution due to GW losses, we can use equation 
([6]) for fc. We must also reverse the sign of (pL because, as 
mentioned above, the GW drag forc e is really applied to 
M2 and not to the test particle as in IMurravf (|l994j ) . We 
now have an expression for the location of the Lagrange 
points in the presence of GW losses: 



sm (pL 



(2-2cos0l)3/2 



(9) 



This relation is plotted as a dashed curve in Figure [21 In- 
terestingly, to leading order, the location of the Lagrange 
points is a function only of the binary orbital separation, 
not the mass ratio, even though the inspiral timescale 
is strongly dependent on /i2. Note that, when including 
non-conservative forces such as GW evolution, there is a 
lack of symmetry between L4 and L5, with L4 moving 
towards the secondary mass, while L5 moves away from 
it. We also find differences in the relative stability of 
orbits around each Lagrange point, as described in the 
next section. 

3. NUMERICAL TESTS 

In Figure [5] we plot the position of test particles near 
L4 and L5 as a function of the binary separation a for 
a range of mass ratios. The initial binary separation is 
ao = 40M. For the most part, the test particles remain 
very close to their initial locations throughout inspiral, 
evolving according to equation © . While our Newtonian 
equations of motion do not formally predict an inner- 
most stable circular orbit (ISCO), we expect the BHs to 
plunge and merge around a = 5M, at which point very 
little evolution in (pL has taken place. The slight offsets in 
the various curves at large a are due to the fact that is 
measured with respect to the origin, not Mi . The inclu- 
sion of post-Newtonian (PN) corrections to the equations 
of motion would likely change the shape of the curves in 
Figure [2] slightly, as well as the detailed shape of sta- 
bility regions around L4 and L5 (iRosswog: & Trautman] 
|1996[ ). It should be noted that Seto & Muto (2010) in- 
clude only IPN terms, and do find a significant reduction 
in stability for certain mass ratios. Yet due to the noto- 
riously uneven converg ence of most PN expansions (e.g. 
iBuonanno et al.l ()2009f ) and references therein), it is not 
clear whether their inclusion would add significant phys- 
ical insight to the problem. 

When studyi ng inert i al dra g forces such as those de- 
scribed above, IMurravl (|1994( ) found that the L4 point 
was more stable than L5. Since our problem is the 
mirror-image of that one, it is not surprising then that 
we find greater stability around L5 when including GW 
evolution. As shown in the upper panel of Figure [2] for 
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Figure 2. Upper: The location of test particles initially at rest 
near the L4 Lagrange point in the rotating frame. As the binary 
system evolves to smaller separation a, the L4 point moves towards 
the secondary BH, located at </) = 0. For fi2 = 0.038, L4 becomes 
linearly unstable around a = 15M, and the test particle is ejected 
at a = lOM. Lower: The location of test particles initially at rest 
near the L5 Lagrange point. Towards the end of the BH inspiral, 
the L5 point moves away from the secondary. In both figures, we 
plot the analytic prediction of equation l(9j as a dashed curve. 



fi2 = 0.038, a test particle initially very close to L4 will 
eventually get ejected from the equilibrium configuration 
as the binary approaches merger. On the other hand, 
test particles librating around L5 appear to become more 
stable, asymptotically app roaching the exac t equilibrium 
position. As mentioned in lErdi et al.l JIQQQ'), for binaries 
with /i2 slightly above the critical value for linear stabil- 
ity (0.0385 < At2 < 0.0401), there still exist bound orbits 
that librate around L4 and L5 with moderate amplitude. 
We find that for these systems, the inclusion of radiation 
reaction increases the stability of L5 by damping the li- 
bration amplitude. On the other hand, the libration am- 
plitude of particles around L4 grows rapidly towards the 
end o f inspiral, as i n the h near case in Figure [5] In con- 
trast, [Seto^&iSut3 (|2Q10f ) find no difference between L4 
and L5, although they only consider a few test particles 
around each Lagrange point. 

To further investigate this asymmetry in stability, we 
carried out orbital integrations for large numbers of test 
particles distributed in clouds around each Lagrange 
point for masses ^2 = 0.01, 0.02, and 0.038. This al- 
lowed us to probe the regions relatively far from the 
equilibrium points, where standard linear stability anal- 
ysis breaks down. We used two methods for estimating 



the large-scale stability of orbits around each Lagrange 
point. The first was a simple number count: the fraction 
of particles that remain bound around the equilibrium 
point. In this context, a particle is considered bound if 
its instantaneous Jacobi constant Cj is less than Cjo(i3) 
for the present orbital separation, corresponding to the 
tadpole-shaped regions around L4 and L5 in Figure [TJ 
Additionally, any particle that escapes with r > 10a or 
approaches closer than 3/zi from the primary or 3/^2 from 
the secondary is removed from the system. The second 
method of measuring stability was to calculate the total 
volume filled by bound particles around each Lagrange 
point. A larger volume corresponds to lower stability, as 
the particles begin to escape the bound libration region. 

The fraction of bound particles /bound is plotted in 
the top panels of Figure [31 normalized to the bound 
number at binary separation of a = 35 M. The parti- 
cles begin in a small cloud of radius SR around each 
Lagrange point with random velocities, isotropic with 
(v) = Voi-hiSR/aY^^ in the corotating frame. We se- 
lect SR to be large enough so that the cloud of parti- 
cles roughly fills the phase-space volume of stable or- 
bits around L4 and L5. In practice, this typically corre- 
sponds to half of the particles getting ejected on roughly 
a dynamical time after the start of the simulation (hence 
/bound is defined relative to a = 35M, not a = 40 Af). Af- 
ter a stable cloud is formed around each Lagrange point, 
marginally bound particles are slowly ejected from the 
system (or accreted onto one of the BHs) as the binary 
system evolves adiabatically. From Figure [21 we might 
expect more particles to get ejected from L4, which ap- 
pears to be less stable, at least at linear order. However, 
as evident in Figure [31 L4 retains more bound particles 
when fi2 = 0.01, L5 appears more stable for fi2 — 0.02, 
and they are nearly identical for /i2 = 0.038. 

It is not entirely clear why this peculiar behavior oc- 
curs, with no apparent trend in global stability as a 
function of mass ratio, nor why one Lagrange point or 
the other should be more stable for a given mass ratio. 
As described in iMurrav &: Dermot3 (jl999[ ). test parti- 
cles near the stable Lagrange points move along epicy- 
cles with two distinct frequencies. For certain values of 
/i2, these libration frequencies have integer commensu- 
rabilities, w hich naturally leads to st rong resonant in- 
teractions fP eprit fc Deprit-Bartholom c 1967). For ex- 
ample, when /i2 = 0.0243, both L4 and L5 are stable 
under standard linear analysis, yet the two epicyclic fre- 
quencies have a ratio of 2 : 1, leading to resonant ex- 
citations and eje ction of all test particles after a few - 
dynamical tim es (|Erdi et al.l 120091 : iSeto fc Mutoll2010( ). 
iMurravl ()1994[ ) shows that these characteristic frequen- 
cies change in the presence of inertial drag forces, which 
we have shown above to be directly analogous to GW 
evolution. Thus we suspect that the different behavior 
seen in /bound (i) for different Lagrange points is due to 
a combination of resonant interactions and asymmetric 
shifts in the libration frequencies due to the binary evo- 
lution via GW losses. 

In the bottom panels of Figure [3] we plot the total 
volume Vbound occupied by bound particles, normalized 
to the overall scale of the binary a^. In the limit of 
a very small 6R/a, where the test particles are well 
within the linear regime, the cloud around each Lagrange 
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Figure 3. Upper: The fraction /bound of test particles that remain bound in Hbrating orbits around L4 {solid curves) and L5 {dashed 
curves) for a range of mass ratios. The particles are initially distributed in a small cloud around each Lagrange point with binary separation 
a = 40M. /bound is normalized to the number of bound particles when a = 35M. Lower: The volume of space occupied by bound test 
particles around L4 {solid curves) and L5 {dashed curves), normalized by the characteristic volume of the system, a^. 



point will actually expand adiabatically relative to the bi- 
nary separation, although shrinking on an absolute scale. 
[Fleming fc Hamilto n ( 2000) show this evolution to scale 
as a^/'': 

where Vi and Vf are the initial and final volumes of test 
particle clouds. We have confirmed this result numer- 
ically with simulations starting with SR <C a. On the 
other hand, if SR/a is large enough such that the stable 
region of phase space is completely filled, the normalized 
volume V/a^ should be constant, since the entire prob- 
lem is scale- invariant. Any adiabatic expansion will not 
increase the volume of bound orbits, but rather decrease 
the total number of particles in these regions, as seen in 
the upper panels of Figure |3l 

We actually find a gradual increase in Vbound/a'^ 
as the binary shrinks, yet still much slower than the 
[Fleming fc Hamilton (2000) resuh of (V/a^) - a'^/'^, 
suggesting that the bound region in phase space is in fact 
filled. The increase in Vbound/a^ during the adiabatic in- 
spiral phase (a > lOM) is likely due to a combination of 
effects. First, we expect a net increase in global stability 
caused by the GW losses, analogous to the effects of in- 
ertial drag described above. Secondly, there is also some 
inherent error in measuring the volume of bound orbits, 
based on our crude metric of the test particles' instanta- 
neous Jacobi constants. Some fraction of particles that 
have moved outside the stable region may still appear 
bound for a few dynamical times before getting ejected 
or accreted, and thus artificially contribute to Vbound- At 
smaller a, where the libration time becomes comparable 
to the inspiral time, the system no longer evolves adi- 
abatically, and a cloud of particles can get "frozen" at 
roughly constant volume, leading to the rapid increase 
in Vijound/a^ for a < lOM. Interestingly, we do see a 



consistent trend of larger Vbound around L4 than for 
all values of /i2, as suggested by Figure [2j However, this 
may be less a measure of local stability than a conse- 
quence of the fact that L4 moves toward the secondary, 
while L5 moves away from it. Thus the two regions be- 
gin to sample different potentials, and the gravitational 
force of the secondary acts to disrupt the cloud around 
L4 more strongly than that around L5, in turn leading 
to a larger Vbound- 

We should add as a word of caution, that many of these 
numerical stability results may be altered somewhat by 
the incl usion of additional P N terms to the equations of 
motion (jSeto fc Mutoll20l"o| ) . In general relativity, unlike 
the case of Newtonian gravity, test particles moving in 
the Kerr metric have three distinct epicyclic frequencies. 
These additional frequencies may lead to more compli- 
cated resonant interactions and thus change the behav- 
ior of /bound(a) and Vbound(a), at least at a quantitative 
level. Yet our intuition suggests that the qualitative be- 
havior shown in Figure [3] should remain unchanged. Fu- 
ture work will investigate this question in greater detail. 

4. ASTROPHYSICAL APPLICATIONS 

4.1. Formation Mechanisms 

In the previous Sections, we showed analytically and 
numerically how the location and stability of the L4 and 
L5 Lagrange points evolve in the presence of radiation 
reaction. We found that test particles, once captured 
into libration orbits around one of the stable Lagrange 
points, will in general remain there throughout the adi- 
abatic inspiral phase, as the entire length scale of the 
system steadily shrinks. However, the initial premise 
was not addressed: will there even be gas or stars or 
debris present in the first place? How might matter get 
captured into these orbits at large binary separations of 
a > 1 pc? Furthermore, are binary BHs with /i2 < 0.04 
astrophysically likely, or even possible? 
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Even in a system as well-studied as the Sun- Jupiter bi- 
nary, there is still no clear consensus on the origin of the 
Trojan asteroids, despite the fact that orbital elements 
are no w known for over 4000 individual objects. iRabd 
([1972) pro posed t hat they may be captured comet s, while 
[Shoema ker et all ()1989f ) and IKarv fc Lissaueri Hl995[ ) 
suggested a more local origin as planetesimals in the early 
solar system. In most theories, some form of dissipativ e 
force is required to capt ure the Trojans (lYodeii I1979D. 
such as nebular gas drag (jMurrav 1994; Karv & Lissauer" 
[l995) or planetesimal collisions (Shoemaker et al. 1989). 
Additionally, an increase in Jupiter's mass via accretion 
of the sur rounding disk could lead to asteroids getting 
captu red (|Marzari fc Scholll [19981 [Fleming fc Hamilton! 
[200l . Analo gous to our BH binary system, the inward 
migration of Jupiter may have also affecte d the capture 
and su bsequent evolution of the Trojans (lYoderi [19791 : 
[Fleming fc Hamihon .2000) . 

The characteristic separation at which point a hard BH 
binary (i.e., the orbital potential is dominated by the 
two BHs, and not by the surrounding stars) is formed 
is comparable to the primary BH's radius of influence 
p/[erritt,2006,l : 



where q — MijMi and 

GM ( M \ f a \-2 

^infl = — - 1-1 PC (^Y^ ) i200kms-J ' 

with a being the one-dimensional velocity dispersion of 
the galaxy's bulge. In the limit of small /i2 ~ g ^ 1, we 
have 

^ VG.02y \ WMq) V200kms-i/ ' 

(13) 

corresponding to ah ~ lO^Af in geometrized units. 
Within this region, t here arc typically thousands of stars 
(jHopman fc Alexande r 2006). and in active galactic nu- 
clei (AGN), a great deal of gas and dust as well. It is 
also roughly the separation at which point the binary 
will merge wi thin a Hubble time due to GW losses alone 
(|Petersiri964[ ). and thus equations ([Ta|) - (lldl) are appro- 
priate. 

For AGN with massive accretion disks, it is 
possible for significant star formation to occur in 
the region where the disk becomes self-gravitating 
(iLevin fc B eloborodov' '200^ Fo r standard thin disk 
models IShakura & Sunvacv (|1973[ ). this self- gravitating 
radius actually corresponds closely to ah. Since galaxy 
mergers are often associated with gas inflow and quasar 
activity, it is quite possible that stars could form in an ac- 
cretion disk just as a hard BH binary is forming, thereby 
trapping the stars around the s table Lagrange points. A 
similar possibility, suggested bv lGoodman fc TanI (|2004l ). 
is that supermassive stars with M > IO^Mq could 
grow out of main-sequence seeds embedded in the ac- 
cretion disk, quickly evolving and then collapsing to 
intermediate- mass black holes (IMBHs) . This scenario is 
particularly attractive, as it naturally provides a mecha- 
nism for making a binary BH with fj,2 < 0.01 on a circular 



orbit, surrounded by ample debris that could easily get 
trapped into resonant orbits. 

Another mechanism that could explain both the low- 
mass secondary and also provide a supply of stars would 
be the tidal disruption of a globular cluster by the 
primary BH. If the glob ular cluster hosts an IMBH 
([Miller fc Hamilton! [2002[ ) . it could form a hard binary 
with the primary, and the surrounding stars could get 
captured into stable orbits around L4 and L5. At this 
point, it would be impossible to carry out any quan- 
titative analysis of the relative rates and likelihoods of 
these different scenarios, considering the remarkable lack 
of observational evidence for even a single SMBH binary 
with orbit smaller than Oh. Rather, we might realisti- 
cally hope to discover such a system with help from one 
of the EM signatures described below. 

4.2. Tidal Disruptions 

Tidal disruptions of solar-type stars by SMBHs are ex- 
pected to occur at a rate of ^ 10~^ yr~^ per MW-type 
galaxy ([Wang fc Mer ritt 2 0(3^ . and are excellent tools 
for detecting otherwise-quiescent nuclei at cos mological 
distances (Lidskii fc Ozernoi!!l979t !Reesl!l988D . Typical 
luminosities are comparable to the Eddington limit, and 
the timescale for the rise and subsequent decay of the 
light curve is on the ord er of weeks to months (jUlmed 
!l999l : [Halpern et al.!!2004l ). which is ideal for detection by 
wide-field time-domain surveys such as Pan-STARR^I 
or LSSTEI ([Gezari et al.!!2009t ). Tidal disruptions have 
recently been proposed as tools for detecting BH merg- 
ers or recoiling mer ger remnants, in some cases long after 
the merger occurs (Kom ossa fc Merritt![2Q08[: [Chen et all 
2009.; Stone fc L oeb 201^^ 

In our scenario, the main-sequence stars are captured 
into resonant motion around L4 or L5 at large separation 
(a ^ ah) and then move in closer to the primary as the 
binary separation slowly shrinks. A star of mass and 
radius will get disrupted when the tidal forces from 
the BH are comparable to the star's own self-gravity, 
which occurs at a distance Rtd from the primary: 

In our numerical tests of Section [Sj we found that during 
inspiral, some fraction of test particles get ejected from 
L4 or and disrupted by M2, the the majority remain 
bound until the final merger stage. For stars that re- 
main bound to L4 or L5, the gravitational tidal forces are 
dominated by the primary, and thus the relevant mass in 
equation ([Hj) is Mi. 

In Figure [1| we plot the orbital separation of a BH bi- 
nary as a function of time before merger {black curves). 
The dashed horizontal lines correspond to the points at 
which (top-to-bottom) an 0/B star, a solar- type star, 
or a white dwarf would get tidally disrupted by the pri- 
mary BH of mass IO'^A/q (top) or IO^Mq {bottom), and 
secondary mass p,2 = 0.01. The vertical dotted lines 
denote where the binary separation reaches the ISCO 
and the BHs merge, and the vertical dashed lines cor- 
respond to the point where the binary enters the LISA 

http: //pan-Starrs . If a.hawall . edU| 
^ http: //mwm. Isst . org| 
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Figure 4. Upper: Tidal disruption radius for a star orbiting a 
10°Mq BH, as a function of time before merger. The secondary 
BH has mass 10^ Mq. Lower: Tidal disruption radius for a star 
orbiting a lO^'M© BH, as a function of time before merger. The sec- 
ondary BH has mass 10^ Mq. FYom top to bottom, the horizontal 
lines correspond to the tidal disruption radius of a main sequence 
O/B star, a solar-type star, and a 0.7 Mq white dwarf. The ver- 
tical dotted line corresponds to the point at which the secondary 
BH roaches the ISCO and the vertical dashed lines the point at 
which the binary enters the LIGO band (/gw ~ 1 mHz). 



GW frequency band. For Mi = IO^Mq, a solar- type 
star would get disrupted roughly a year before merger. 
At that point, LISA would already be able to provide an 
error-box o f ^ 100 s quare degrees for the source location 
on the sky (iLang fc Hughes 2008), which could be easily 
monitored for lumi nous tidal disruptio n flares out to red- 
shifts z ^ 0.2 or so (jGezari eEiDIlOOi). In ah cases with 
Ml > IO^Mq, white dwarfs would not get disrupted, but 
rather plunge directly into the BH. 

In addition to the stars that get disrupted at late times, 
we saw in the simulations of Section [3] that there are also 
some particles that get captured by either Mi or M2 dur- 
ing the inspiral phase, thus leading to tidal disruption, 
or get ejected entirely from the system. Table [1] shows 
the fraction of test particles that meet one of these var- 
ious fates during inspiral, i.e., when the binary separa- 
tion is 5M < a{t) < 35M. Of course, if equation is 
met, then even stars on stable orbits are disrupted. As 
in Figure |31 there are no clear trends with varying ^2 , 
but we do confirm that particles around L5 are generally 
more likely to get ejected from the system, while particles 
around L4 are more likely to get captured or disrupted 



Table 1 

Fraction of test particles that, during inspiral, are either ejected 
from the system or captured by one of the BHs. 



At2 = 0.01 tM2 = 0.02 ^l2 = 0.038 
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by the BHs. 

4.3. Ultra-Hyper- Velocity Stars 

Of the stars that do get ejected from the binary system 
during inspiral, we find a wide range in their velocities at 
infinity. To first order, one would expect the ejection ve- 
locity to be comparable to the orbital velocity. However, 
particles ejected at late times when the orbital velocity 
is highest are also coming from the deepest potentials, so 
are slowed more as they escape from the system. In prac- 
tice, we find that net effect is a slight increase in the av- 
erage escape velocity at late times, but somewhat flatter 
than the orbital velocity scaling of a~^/^. This is likely 
due to the fact that the most efficient three-body sling- 
shot interactions require the test particle to come within 
a small fraction of a from the secondary, and when a be- 
comes smaller, these close approaches more likely lead to 
capture or disruption than ejection. 

From equation (|14p , we see that a solar- type star would 
get tidally disrupted by a IO^Mq BH at a distance of 
~ 50M, far outside the ISCO. For Mi = lO^M©, how- 
ever, a star could survive to within lOAf of the primary. 
From that separation, we find that test particles are 
ejected with mean velocities of 0.1c, roughly 50 times 
greater than the radial velocities of hyper-velocity stars 
HVSs) observed in the outer regions of the Milky Way 
Brown et al.ll2009h . Of course, any star moving that 
fast would be observable for a much shorter time; only 
a few million years after the BH merger. At the same 
time, this allows younger, more massive stars to survive 
out to the edge of the galaxy, making them even easier 
to identify against the field stars. It is conceivable that 
the BH at the center of our galaxy may have experienced 
a minor merger with fi2 ^ 0.01 in the past 10 Myr, an 
event which may be convincingly inferred by the detec- 
tion of such an ultra- HVS at a distance of ^ 50—100 
kpc. 

4.4. Compressed Stellar Systems 

In the event that many stars are initially trapped 
around L4 or L5, the gravitational interactions between 
the test particles may become important, an effect that 
we have neglected in this work. Especially for larger sys- 
tems with M > IO^Mq, one could even imagine forming 
globular clusters around the Lagrange points, although 
from equation it is clear they would necessarily 

be quite compact. When the core of a normal globu- 
lar cluster grows too dense, it can undergo a collapse 
via gravitational instability, at the same time expanding 
the orbits of the outer stars, due to conservation of en- 
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ergy. In our BH binary system, the adiabatic compres- 
sion of orbits may trigger core collapse of the trapped 
cluster, while simultaneously preventing any outer stars 
from escaping, thereby enhancing the stellar density fur- 
ther. If the density becomes high enough, stellar col- 
lisions will become common, and the system could act 
as a breeding ground for sup ermassive stars and IMBHs 
(jPortegies Zwart et al.|[l999[) . 

4.5. Compact Objects 

Many of the formation scenarios considered above nat- 
urally lead to high-mass stars forming in the dense 
gaseous environments around the SMBH binary. These 
stars could then evolve and collapse into neutron stars 
or stellar-mass BHs on timescales short compared to the 
inspiral time. Unlike main sequence stars that would get 
tidally disrupted well before the merger, compact objects 
would survive unperturbed until they either get ejected 
or captured whole by one of the BHs. In that case, they 
would also contribute to the GW signal produced by the 
binary. For a system with Mi = lO^M©, M2 = IO^Mq, 
and mo = lOM©, the stellar- mass BH would contribute 
an additional 0.1% to the waveform amplitude. At a dis- 
tance of 1 Gpc, this small perturbation would still be 
observable with LISA. 

However, if the compact object is orbiting exactly at 
the L4 or L5 with no libration, it would contribute a 
pure sinusoid component to the waveform with identi- 
cal frequency as the dominant wave, only with a phase 
shift of ±120°. Thus it would be indistinguishable from 
a normal binary on a circular orbit with slightly different 
masses and phase. On the other hand, if the compact ob- 
ject is librating with typical amplitude w 30°, as we 
might expect in the general case, the waveform will con- 
tain a distinct signal at the libration frequency (in this 
example, wnb ~ 0.27cJorb) which may still be detectable. 
More likely, this extra signal would be confused with a 
spin-orbit or spin-spin precession frequency, effectively 
degrading the overall confidence in the binary parame- 
ter estimation. In some cases, the stellar-mass BH might 
get ejected during inspiral, but with small enough ve- 
locity that it doesn't escape the system, but rather re- 
mains bound on a larger, eccentric orbit even after the 
two massive BHs merge, ultimately merging much later 
as an extreme mass-ratio inspiral (EMRI) source. 

As mentioned above, the secondary BH and the "test 
particle" can actually have comparable masses and still 
maintain stability of the Lagrange points, as long as 
/X2 < 0.0385. Thus a stable three-body BH system could 
exist with Mi w IO^Mq and Ah « M3 « lO'' Mq, in 
which case the orbital libration could be seen very clearly 
in a LISA waveform. The probability of such a system 
existing is perhaps not as small as one might think, espe- 
cially if both IMBHs form out of a massive accretion disk, 
not unlike Saturn's moons Janus and Epimetheus. In- 
deed, only 15 years ago, the possibility of "Hot Jupiters" 
with /i2 > 0.001, orbiting very close to their host stars 
was considered quite remote, and now we know of hun- 
dreds of such systems. The same may very well be true 
of SMBHs and their surrounding gas disks. 

4.6. Diffuse Gas 

In the circumbinary accretion disk model, ample gas 
and dust will surround the BH binary, of which some 



fraction could likely get trapped into orbits around L4 
and L5. If cold enough, this gas could get compressed 
by the evolving binary potential and form stars or even 
clusters of stars, as discussed above. More likely, the 
gas will be rather hot, with thermal velocities compara- 
ble to the libration velocities of test particles filling the 
stable volume around L4 and L5. One very rough esti- 
mate of this thermal motion may come from inspection 
of equation ([2]) and Figure [TJ Trapped particles will have 
Cjo(-^4) ^ Cj < Cjo(i3), so typical velocities in the ro- 
tating frame of 

^ CjiL3,v ^ 0) - Cj{U,v ^ 0) ^ 2^i2a-^ (15) 



( Murrav &: DermottI Il99i^ . Combining with equation 
( T3|) for the separation at the time of binary formation, 
we get an initial thermal velocity of 



Vthc 



600 kms" 



If ^ ) 

V200kms-i/ 



(16) 



or roughly 4 x 10^ K. Somewhat remarkably, Wthcrm is 
independent of fi2 or AIi (except indirectly through the 
M — a relationship). 

We can estimate an upper limit on the density of this 
gas po by requiring it to be stable to gravitational col- 
lapse. In other words, the total mass should not exceed 
the Jeans mass: 



0.005 poafi < 



5fci3T 
Grrir, 



3/2 



1/2 



(17) 



with fcs the Boltzmann constant and rup the proton 
mass. The leading factor of 0.005 comes from the vol- 
ume of the stable region, estimated from Figure [3l Tak- 
ing mpWt\gj„ = iksT, we get 



Po < 35 



"the 



Gal 



(18) 



corresponding to a total mass in gas of 



O.OObpoal < 7 X IO^Mq (^) ( 



M 



lO^MfT 



(19) 

which is far more than enough material to fuel a bright 
EM counterpart. 

Even if the density of the trapped gas is orders of mag- 
nitude smaller, it should still be sufficient to produce ob- 
servable emission lines if irradiated by an AGN or quasar 
powered by accretion onto one of the BHs. For a cloud of 
gas filling the stable region around each Lagrange point, 
the covering fraction would be on the order of a few per- 
cent. Since the thermal velocity is roughly a factor of 
10 smaller than the orbital velocity, the emission lines 
would be somewhat narrower than a classical broad line 
originating from the same radius. The coherent shape 
and motion of the cloud would also give rise to a large 
velocity offset relative to the host galaxy or narrow emis- 
sion line clouds residing at a greater distance from the 
central engine. This offset in wavelength would vary pe- 
riodically with the binary orbit (Torb at our fiducial value 
of afi would be on the order of 10 years), giving strong 
evidence for the existence of a BH binary. Such an ef- 
fect could be detected by large spectroscopic surveys like 
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the Sloan Digital Sky Survey, a technique already in fre- 
quent use by observ ers searching for bina r y quasars or re- 
coilin g black holes (iBonninget alJ [20071 : iKomossa et alj 
[200l lBoroson fc Laueril200l 

5. DISCUSSION 

We have analyzed the circular restricted three-body 
problem including leading-order terms for gravitational 
wave losses, in particular studying the behavior of test 
particle orbits near the stable Lagrange points L4 and 
L5. We find that these orbit remain stable throughout 
most of the adiabatic inspiral phase, and their locations 
shift relative to the Newtonian limit as the binary en- 
ters the late phases of inspiral and merger. Analogous to 
the case of inertial drag forces in the planetary problem, 
the non-conservative GW loss terms lead to an asymme- 
try between L4 and L5: test particles around L4 move 
closer to the secondary BH, while those initially around 
L5 move away from it. While there appear to be compli- 
cating resonances that affect the relative stability of the 
two regions differently for different mass ratios, a general 
trend is that the inclusion of GW losses increases the sta- 
bility of L5 while decreasing the stability of L4. Particles 
initially around L4 are more likely than those around 
to get captured or disrupted by M2, while particles near 
L5 are more likely to get ejected from the system entirely. 
However, in both cases, the vast majority of particles get 
neither ejected nor captured until the binary separation 
is on the order a < 5Af, at which point the Newtonian 
equations of motion break down entirely and everything 
plunges into the primary BH. 

The stability of the L4 and L5 Lagrange points 
throughout most of the inspiral provides a natural mech- 
anism for shepharding material from relatively large radii 
down to the central BH, in turn producing a luminous 
EM counterpart to the GW signal from a BH merger. 
The detailed physics behind such an EM signal are 
well beyond the scope of this paper, but early results 
from numerical relativity simulations suggest that bright 
EM flares sho uld be quite general, provided enough 
gas is a vailable ([ van Meter et al.l[20TotlBode et al.l[2010l : 
iFarris et al.. ,2010l ). In addition to the prompt emission 
associated with the final merger, we have shown that 
hyper- velocity stars with speeds as high as 0.1c might 
be observable, flying away from the galactic center for 
millions of years after the merger. With sufficient col- 
umn densities of diffuse gas trapped around ^4 and L5, 
a precursor signal could take the form of highly red/blue- 
shifted emission lines, narrower than typical broad lines, 
with the line offset varying periodically with the BH or- 
bit. 

For the majority of the results in this paper, we require 
that the secondary mass be relatively small: /X2 < 0.0385 
to guarantee linear stability of L4 and L5. It is unclear 
how common BH binaries with such mass ratios might 
be, and whether they are more likely to form via galac- 
tic mergers or in situ, via the collapse of supermassive 
stars, accretion from a surrounding disk, or capture from 
a tidally stripped globular cluster. However, there are 
other resonances in the restricted three-body problem 
that may lead to qualitatively similar behavior: the sta- 
ble transport of matter from large radii down to the ISCO 
at the time of merger. One example from the solar sys- 
tem is the Hilda family of asteroids that occupy stable 



eccentric orbits with periods 2/3 that of Jupiter. This 
ensures that they avoid close encounters with Jupiter 
at apocenter, where they spend most of their time. In 
the corotating frame, a snapshot of their positions will 
appear to fill out a triangle with vertices align ed with 
L^_L4, and L5 (Broz & VokrouhhckY 2008.). iGueronI 
((2OO6I) showed that this resonance is stable even in the 
extreme relativistic limit, but still with small mass ratios 
(/i2 ~ 0.001). Our own preliminary results suggest that 
this 2 : 3 resonance is linearly stable for mass ratios at 
least as large as /i2 ~ 0.1, but we leave a comprehensive 
study of the problem, including GW losses, to a future 
paper. 

In addition to all the Newtonian resonance behavior 
present in the three-body problem, the two BHs may also 
get locked into spin-orbital resona nces with each other 
during inspiral (jSchnitt man' '2004). Coupled with the 
non-degenerate orbital frequencies of test particles in a 
Kerr background, the inclusion of spinning BHs would 
introduce many new degrees of freedom that may affect 
the stability of L4 and L5. While we do not think ad- 
ditional PN terms in the equations of motion will signif- 
icantly change our results on a qualitative level, future 
work that addresses these questions in more detail is cer- 
tainly warranted. 

The author would like to thank Doug Hamilton, Matthew 
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work was supported by the Chandra Postdoctoral Fel- 
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